/*******************************************************************************
 * Copyright 2009 DCSpectrometer - http://code.google.com/p/dcspectrometer 
 *  
 * Licensed under the Apache License, Version 2.0 (the "License"); 
 * you may not use this file except in compliance with the License. 
 * You may obtain a copy of the License at 
 *  
 *     http://www.apache.org/licenses/LICENSE-2.0 
 *     
 * Unless required by applicable law or agreed to in writing, software 
 * distributed under the License is distributed on an "AS IS" BASIS, 
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 
 * See the License for the specific language governing permissions and 
 * limitations under the License. 
 *******************************************************************************/
package com.dcspectrometer.Analysis;

/**
 * This class calculates Cumulative Distribution Function for Normally Distributed Variables
 * Class is static and cannot be instantiated. 
 */
public final class NormCDF {
	private static final double coefficient = 2 / Math.sqrt(Math.PI);

	/**
	 * Calculates probability that normal random variable is less than x
	 * @param x - input value
	 * @return - Probability that normal random variable is less than x
	 */	
	public static double getCDF(double x) {
		double tolerance = 1e-15;
		int maxCycles = 100;
		return 0.5 * (1.0 + getERF(x/Math.sqrt(2.0), tolerance, maxCycles));
	}
	
	/**
	 * Estimates error function value for x
	 * @param x - input value
	 * @return - error function value for x
	 */
	private static double getERF(double x, double tolerance, int maxCycles) {
		int n = 1;
		double nextCorrection = tolerance + 1;
		double erf = x;
		while((n <= maxCycles) && (Math.abs(nextCorrection) > tolerance)) {
			nextCorrection = Math.pow(x, 2*n + 1)/(factorial(n)*(2*n + 1));
			if((n % 2) == 1) erf -= nextCorrection;
			else erf += nextCorrection;
			n++;
		}
		erf = erf * coefficient;
		return erf;
	}
	
	/**
	 * Calculates n factorial.
	 * @param n - input value
	 * @return
	 */
	private static double factorial(int n) {
		double result = 1.0;
		for(int i=2; i<=n; i++) {
			result *= i;
		}
		return result;
	}
}
